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We report the discovery and confirmation of a transiting circumbinary planet (PHI) 
around KIC 4862625, an eclipsing binary in the Kepler field. The planet was discovered 
by volunteers searching the first six Quarters of publicly available Kepler data as part 
of the Planet Hunters citizen science project. Transits of the planet across the larger 
and brighter of the eclipsing stars are detectable by visual inspection every ~137 days, 
with seven transits identified in Quarters 1-11. The physical and orbital parameters 
of both the host stars and planet were obtained via a photometric-dynamical model, 
simultaneously fitting both the measured radial velocities and the Kepler light curve of 
KIC 4862625. The 6.18 it 0.17 planet orbits outside the 20-day orbit of an eclipsing 
binary consisting of an F dwarf ( 1.734 ± 0.044 R©, 1.528 ± 0.087 Mq) and M dwarf 
( 0.378 ± 0.023 Rq, 0.408 ± 0.024 Mq). For the planet, we find an upper mass limit 
of 169 M© (0.531 Jupiter masses) at the 99.7% confidence level. With a radius and 
mass less than that of Jupiter, PHI is well within the planetary regime. Outside the 
planet's orbit, at ~1000 AU, a previously unknown visual binary has been identified 
that is bound to the planetary system, making this the first known case of a quadruple 
star system with a transiting planet. 

Subject headings: binaries-eclipsing — stars: individual (KIC 4862625) — planets and satel- 
lites: detection — planets and satellites:general 



1. Introduction 



For the past three and a half years NASA's Kepler spacecraft (iBorucki et al.ll2010l ) has mon- 
itored, nearly continuously, over 160,000 stars looking for the signat ures of transiting exo planets. 
The Kepler team has discovered more than 2300 planet candidates (jBatalha et al.ll2012l ). nearly 
quadr upling the samp l e of known planets, in addition to identifying close to 2200 eclipsing binary 



stars ( Prsa et al.ll201ll : ISlawson et al.ll2011 



One of Kepler's many successes has been the first con- 
firmed i dentification of transiting circumbina r y planets, planets o rbiting both stars in a binary star 
system (jPovle et al.l boill : IWelsh et al.l I2OI2I : lorosz et al.l l2012al lbh . Although plan ets have previ- 



ously been suspected to be orbiting both stars in a binary star s ystem before Kepler (jBennett et al 



1999 



Deeg et allboosl : iLee et~al]l2009l : lieuermann et al.lboid v the Kepler era ushered in the first 



confirmed identification of such systems. 

Identifying transiting exoplanets in binary star systems, on P-type (circumbinary) orbits where 
the planet orbits around both stars, is difficult when using the automated transit search methods 
and ca ndidate validation processes employed by the Kepl er team. The Transiting Planet Search 
(TPS) (jTenenbaum et al.ll2012l : Ijenkins et al.ll2010bl . I2OO2I ) pipeline and the c andidate vetting pro - 
cess were primarily designed to identify transiting planets around single stars (jBatalha et al.ll2012l ). 
Transits due to planets in binary star systems can easily be drowned out by the much larger stellar 
eclipse signal and may even be superimposed on top of a primary or secondary eclipse escaping 
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detection. What makes automatic detection challenging is the fact that the planet transits in such 
systems are not necessarily constant in duration or repeated at regular intervals due to the chang- 
ing positions and velocities of the stars at each transit. Alternative detection techniques, including 
visual inspection and the search for eclipse timing variations on the host stars due to non-Keplerian 
three-body effects, are currently being employed to find planets in the nearly 750 eclipsing binary 
systems with orbital periods larger than 1 day. These efforts have proven successful with the discov- 
ery of six transiting circumbinary pl anets in the Kepler field (Kepler-16b , Kepler-34b, Keple r-35b, 
Kepler-38b, Kepler-47a, KeDler-47b: loovle et al.ll201ll : Iwelsh et allboii lOrosz et al.ll2012al lbh. 



These multistar planetary systems are the extremes of planet formation and serve as the test 
grounds for planetary accretion and evolution models th at must be able to explain the e xistence of 
planet s in these dynamically challenging environments. iPaardekooper et al.l (|2012l ) and iMeschiari 
([2013) are unable to simulate the creation of planetary cores at the distances the Kepler circumbi- 
nary planets currently reside suggesting these planets formed further from their parent stars and 
migrated inward to their present orbits. Models of circumbinary protoplanetary disk evolution by 
Alexander! (j2012l ) predict that disks around close binaries with semimajor axes (a) less than 1 AU 
live longer than those around single stars, but disk lifetimes decline as photoevaporation increases 
at larger binary semimajor axes. As a result, they predict a dearth of circumbinary planets around 
wide b inaries with a > 1 AU and an abundance of circumbinary planets in stellar binaries with a < 
1 AU. iBate et al.l (j200d ) predict that circumstellar disks around binaries with a < 100 AU should 
be coplanar to the orbit of the binary. Thus any planets formed from such a disk are expected to all 
be ali gned with the binary's orbit. According to dynamical evolution models by lPierens &: Nelson 
(|2008l ). Jupiter-mass planets on circumbinary orbits will be rare and difficult to find. In their sim- 
ulations, the Jupiter-mass bodies typically had a close encounter with the secondary star during 
migration causing the planet to be subsequently ejected from the system or scattered outward to 
the edge of the circumstellar disk. Therefore planets with masses equal to or less than that of 
Saturn are predicted to be the most abundant in these types of planetary systems. By finding more 
planets in multistar systems we can probe the environments these planets reside in and begin to 
test and evaluate the predictions from these formation and dynamical evolution models. 

We report here the discovery and characterization of a new transiting circumbinary planet, 
orbiting a known Kepler eclipsing binary, via visual inspection of the first 16 months of publicly 
released Kepler data by volunteers from the Planet Hunters citizen science project. This is Planet 
Hunters' first confirmed planet, and we subsequently refer to the planet as PHI'. The 6.18 =b 0.17 

planet orbits outside a 20-day period eclipsing binary. Beyond the orbit of PHI, a distant 
visual binary located at a projected distance of ~1000 AU is bound to the planetary system. This 
is the first identification of a confirmed transiting planet orbiting in a quadruple star system. In 
this Paper, we present the discovery, observational follow-up, and characterization of this hierar- 
chical stellar-planetary system. Combining the Kepler photometric light curve and radial velocity 
observations with a 3-body photometric-dynamical model, we constrain the physical and orbital 
properties of the planet and host stars. 
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2. Discovery 



The Planet Hunter^ ( Fischer et al. 2012 : Schwamb et al. 2012 ) citizen science project uses 
the power of human pattern recognition via the World Wide Web to identify tr ansits in the Kepler 



public data. Plan et Hunters uses the ZooniverscI platform, first described in ( Lintott et al. 



2008 



Smith et al.ll201ll ). to present visitors to the Planet Hunters website with a randomly selected ~30- 
day light curve segment from one of Kepler's ~160,000 target stars. 5-10 independent volunteers 
review each 30-day light curve segment and draw boxes directly on the web interface to mark the 
locations of visible transits. For further details of the website and classification interface we refer 



the reader to (ISchwamb et al.ll2012l ). 



In addition to the main classification interface, the Planet Hunters website hosts an object- 
oriented discussion and investigation tool known as 'Planet Hunters Talk' (referred to as 'Talk'jfj. 
Talk hosts forum-style message boards enabling volunteers to discuss and analyze interesting light 
curves together. In addition each 30-day light curve presented on the Planet Hunters main classifi- 
cation interface has a dedicated page on the Talk website where volunteers can write comments, add 
searchable Twitter-like hash tags, and group similar light curves together. Volunteers are directed 
to these pages through the main Planet Hunters website. Once a classification on the main interface 
is complete, the volunteer is asked "Would you like to discuss this star?" . If the volunteer answers 
'yes', he or she is directed to the light curve's dedicated Talk page. From there a Planet Hunters 
volunteer can also view all available public light curve data for that Kepler target presented as it 
would be seen on the main Planet Hunters website, including the ability to zoom in both time and 
brightness. 

Talk was designed to aid in identifying objects of particular interest and unusual light curves 
that are difficult to find via automatic detection algorithms and beyond the scope of the main 
Planet Hunters classification interface. It also provides an avenue for communicating directly with 
the Planet Hunters science team through the message boards and comments. From 2010 December 
16 to 2012 June 9, in total, 69 percent of registered Planet Hunters participants have visited the 
Talk site with 23 percent having contributed comments. In this time frame, an average of 116 and 
a median of 72 Planet Hunters users visit Talk each day. In addition an average of 511 and median 
of 329 comments discussing presen ted light curv e s are generated on Talk daily. We note that the 
two planet candidates presented in lLintott et al.l ( 2012 ) were first identified by the Planet Hunters 
science team via the discussion of the light curves on the Talk message boards before a systematic 
search of the Q2 light curves had been completed. 



On 2011 May 11, coauthor Kian Jek posted on Talk highlighting the light curve of eclipsing 



^http: / /www.planethunters.org 
^http: / /www.zooniverse.org 

^http://talk. planethunters.org - The code is available under an open-source license at 
https://github.com/zooniverse/Talk 
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binary KIC 12644769, identifying additional transits due to a possible third body. ISlawson et al, 



(|201ll ) had already noted the additional transit features, and the sys tem subsequently turned out 



to be Kepler-16, the first circumbinary discovered in the Kepler field (|Dovle et al.ll201ll ). After the 



announced di scovery of Kepl e r- 16b , the Planet Hunters science team posted links to the Talk pages 
for the entire ISlawson et alj (j201ll ) eclipsing binary list encouraging volunteers to examine these 
eclipsing binary light curves and perform their own search for additional transits due to an orbiting 
planet. During visual inspection of the light curve for detached eclipsing binary KIC 4862625, 
coauthor Robert Gagliano spotted two transit-like features in Quarters 2 and 4 (see Figure [T]) 
separated by ~137 days, and on 2012 March 2 highlighted this light curve on the Talk forums as a 
potential circumbinary. At the time only Quarters 1-4 (~310 days of science data) were available 
on Planet Hunters website, but Quarters 1-6 were publicly available. On 2012 March 3, Kian Jek, 
examined the Quarter 5 and 6 light curves, obtained from the Kepler Mission Archivqj hosted by 
the Mikulski Archive for Space Telescopes (mastI), identifying a third transit with similar depth 
in Quarter 5 occurring roughly 137 days after the Quarter 4 transit. These additional transit events 
in Quarters 1-5 as well as others subsequent Quarters were also independently identified by Kostov 
et al.(2012) in a separate and concurrent analysis of the Kepler public data. After examining the 
reported column pixel centroids in the FITS files for a shift during in and out of transit that may 
indicate contamination from a faint background eclipsing binary, Kian Jek noted the presence of 
this third transit in a Talk discussion threajl and notmed the Planet Hunters science team of the 
potential discovery of a circumbinary planet. 



3. Kepler Light Curve 



There are two types of Kepler light curves provided in MAST: Pre-Search Data Conditioning 
(PDC) and Simple Aperture Photometry (SAP). The PDC correction optimizes the light curve for 
detection for planet transits as well as removing systematic effects and long-term trends, but also 
makes some assumption for stellar variability and attempts to remove that signal. The SAP light 
curve is the summed flux in the target aperture with no further correcti ons. For furth e r speci fics 
of the Kepler data proce ssing and data products we refer the reader to Jenkins et al. ( 2010bl lal): 
Christiansen et al.l (j201ll ). and ISmith et al.l (|2012l ). We note that the original discovery on the 
Planet Hunters website was using the PDC-corrected light curve, but for the remainder of the 
analysis presented here, we choose to use the SAP light curves in order to preserve the stellar 
variability and the eclipse signals, performing our own detrending and removal of instrumental 
systematics for each individual Quarter of Kepler observations. 



^http: / /archive. stsci.edu/kepler 
^http: / /arciiive. stsci.edu 



6 



http://talk.planethunters.org/discussions/DPH1014m5m?group_id=SPH10052872 
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The normalized light cur ve for KIC 486262 5, for Quarters 1-11 is shown in Figure [2l The 
Kepler Input CatalogJ (KIC. iBrown et al.l l201ll ). reports a magnitude of 13.72 for the eclipsing 
binary in the broad Kepler filter (Kp). We will refer to the brighter more massive star as 'Aa' and 
use 'Ab' to denote the fainter lower-mass stellar companion in the eclipsing binary. These two stars 
have an orb ital period of 20.0002 14 days on an eccentric orbit with an eccentricity of 0.24181 as 
reported by ISlawson et al.l (120111 ) . In the Kepler light curve the depth of the primary eclipses are 
~1.3%, and the depth of the secondary eclipses are only slightly deeper than the planet transit at 
~0.11%. Three transits of PHI crossing Aa are visible in Quarters 2, 4, and 5, the data publicly 
available at the time of discovery, with another four identified in Quarters 7, 8,10, and 11. No 
transits across Ab were identified. Initial modeling of the light curve to remove brightness changes 
due to stellar variability found that the three transits in Quarters 2,4, and 5 had similar depths 
(~0.1%) consistent with the third body transiting being plane tary in size, ~2.8 Rm assuming a 



primary (Aa) stellar radius of 0.804 Rq as reported in the KIC ([Brown et al 



2011 



One of the largest sources of false positives for Kepler planet candidates are blended faint 
background eclipsing binaries, but this would require a chance sky alignment of two eclipsing binary 
systems along the line of sight. This scenario, though highly unlikely, cannot be completely ruled 
out by statistical arguments alone. If the transiting body truly is orbiting both stars in the binary. 



and velocities of the stars at each transit (e.g lDeeg et al. 



2012 



Deee et al. 


1998; 


Dovle et al. 


2000, 


2011; 


Welsh et al. 



Orosz et al.ll2012al jbh. A transit function is fit individually to each of the seven transits to 



measure the transit times and widths. The measured transit midpoints are reported in Table [H 
The mean period is 135.648 days. We note all times and ephemerides in our analysis are using 
Barycentric Julian Dates in Barycentric Dynamical Time, which we will refer to as BJD. Observed 
(O) minus Calculated (C) offsets and durations for the seven transits are plotted in Figure [H 0-C 
amplitudes as high as it 1 day are observed and the transit durations are varying with the binary 
orbital phase, confirming this body is indeed orbiting stars Aa and Ab and the transit features are 
not due to a blended background eclipsing binary or any other astrophysical false positive. 

Furthermore, the tidal field of the planet may cause nonlinear timing variations in the eclipsing 
binary, but the zeroth order effect is a period difference between primary and s econdary eclipses . 



as had beeri obser ved the for Kepler-16b, Kepler-34b, and Kepler-35b systems (jDoyle et al.l 12011 



Welsh et al.ll2012l ). Primary and second ary eclipse timings were estimated by fitting a low-order 



cubic Hermite polynomial (for details see lSteffen et al.ll201ll : IWelsh et al.ll2012l ) to each individual 
eclipse in the normalized Quarters 1-11 light curve. The measured primary and secondary eclipse 
midpoints for KIC 4862625 are reported in Table [2j The best fitting linear primary and secondary 
eclipse ephemerides are: Tq = -52.18091 ± 0.00018 (BJD - 2455,000) with a period P = 20.0002498 
± 0.0000064 days for the primary eclipse and Tq = -44.3252 ± 0.0016 (BJD - 2455,000) with 
P = 20.000244 lb 0.000056 days for the secondary eclipse. The difference between the primary 



^http: //archive. stsci.edu/Kepler/Kepler fov/search.php 
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and secondary periods is 0.50 it 4.87 seconds. The timings of the stellar eclipses are measured 
and compared to the that expected from the best-fit eclipse ephemerides. The 0-C offsets for the 
primary and secondary stellar eclipses are plotted in Figure [3j No discernible timing variation signal 
in either the primary or secondary eclipses were found. The planet may not be massive enough 
to gravitationally perturb the stellar orbits significantl y over the timescal e of the observations, as 
is the case for the Kepler-38 and Kepler-47 systems ( Orosz et al. 2012al lbl). Without detectable 



perturbations on the stellar orbits, we will only be able to place an upper limit on the mass of PHI 
when the light curve and radial velocity observations are combined in the photometric-dynamical 
model fit as described in Section [71 



Flux Contamination 



The Kepler pix el scale is 4" per pixel, with a typical photometric aperture radius of 6" 
(j Jenkins et al.ll2010al ). With such a wide aperture, stellar contamination and photometric blends 
are a concern when analyzing the Kepler photometry of target stars. Adding linearly, the contribu- 
tion of extra light decreases the observed transit and eclipse depths. Accurately estimating the size 
of the transiting planet requires knowledge of the additional flux contributors to the Kepler light 
curve for KIC 4862625. A search of the KIC found that two nearest so urces are KIC 4862 611 and 
KIC 4862601 located 8.68" and 11.8" respectively from KIC 4862625 (|Brown et ahlboilh . With 
Kepler magnitudes of 18.052 and 18.256 respectively, these two stars would not be completely in 
the Kepler photometric aperture and at most would contribute less than a percent of total flux in 
the light curve. 

Faint companions closer than a few arcseconds to the Kepler ta rget may not be found in 



the KIC or 2MASS catalog (ISkrutskie et al.l l2006l : Brown et al.ll201ll ). To better constrain the 



presence of additional stars not listed in the KIC that may be contaminating the Kepler aperture, 
we obtained near- infrared adaptive optics (AO) observations us ing the NIRC2 im ager behind the 
natural guide-star (NGS) AO system on the Keck II telescope (jchez et ahlboosl l. K'(1.948-2.299 
/im) imaging was taken on 2012 May 29 UT with NIRC2 in the narrow camera setting, yielding 
a 10x10" field-of-view with a pixel scale of 9.9 mas per pixel. We used a 9-point dither pattern 
with a 2-" offset between each position. At each nod, 60s integrations were taken, resulting in a 
total integration time on target of 540s. The raw frames were dark and background subtracted, 
flat fielded. Bad and hot pixels were interpolated over. The reduced frames were then combined 
into a single image. Visual inspection revealed two nearby stars: a source 0.702" southwest of KIC 
4862625 approximately 5 times fainter than the target star and a much fainter source ~3" away 
to the southwest. From this point forward, we refer to these sources as the 0.7" contaminator 
and the 3" contaminator respectively. Both stars are well within the photometric aperture of KIC 
4862625. The Kepler data analysis pipeline would not account for the flux contribution from these 
contaminating stars. 



In order to estimate the light curve dilution due to the brighter 0.7" contaminator, we obtained 
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further observations with NIRC2 NGS AO in Ks (1.991-2.302 /im) and J (1.166-1.330 /xm) filters on 
2012 June 7 UT in the same observing setup. The images were reduced in the same manner as on 
the first night but with a 3" offset between nods. J images were taken with an integration time of 
90s exposures at each nod. The Ks imaging consisted of 60s at each nod. A total integration time 
of 810s in J and 540s in Ks respectively was achieved. There were intermittent high cirrus clouds 
during those observations, and therefore we were only able to obtain relative J and Ks photometry. 
Figure [5] presents the reduced Ks AO image for the 0.7" contaminator and the primary KIC target. 
Using aperture photometry, we measure SKg = 1.67zb 0.03 and 5J = 1.89 it 0.04 mag difference. 



We adopt the method developed by iHowell et al.l (j2012l ) to extrapolate to the broad band 



Kepler filter (Kp) using J and Ks colors. The measured relative J and Ks colors were combined with 
the reported 2MASS photometric magnitudes (J=12.714, K=12.394) to obtain the true apparent 
Ks and J magnitudes for both KIC 4862625 and the 0.7" contaminator. The 2MASS colors are 
a summation of all three identified sources. The 3" contaminator is too faint in AO observations 
to measure an accurate relative magnitude, but contributes at most at the percent level to the 
infrared colors. Therefore, we chose to assume KIC 4862625 and the 0.7" contaminator are the 
only contributions to the reported 2MASS colors . The KIC r eport s Kp= 13.71 mag for KIC 



4862625. Using the filter transforms determined bv I Howell et al.l (|2012l ) we find a Kp magnitude of 
13.79 for KIC 4862625 and 16.10 for the 0.7" contaminator is the Kepler target is a dwarf star. If it 
is a giant, we estimate a Kp magnitude of 13.77 for KIC 4862625 and 16.21 for the faint companion. 
Thus the 0.7" contaminator contributes 11-12% of the light in the light curve for KIC 4862625. 



It is not surprising that we find a source less than 1" from KIC 4862625. I Adams et al.l (120121 ) 



find that 20% of the Kepler candida tes observed with adaptive optics have a close companion 



within 2", and iLillo-Box et al.l (|2012l )report similar findings with 17% of their sample of Kepler 
candidates examined had at least one contaminating source within 3" of the Kepler target. Using 
the KIC stellar parameters and the expected galactic extinction, KIC 4862625 is approximately 1 
kpc away, thus the projected distance between the binary and the 0.7" contaminator is ~1000 AU. 
Although gravitationally bound multistellar systems exist with such separations, we would need 
measurements of the proper motions of both the binary (Aa and Ab) and the 0.7" contaminator 
to confirm association. From the NGS AO observations alone, we cannot determine if the 0.7" 
contaminator is bound to the system. The time separation of our observations is insufficient to 
measure proper motions, but the radial velocity observations (discussed in Section [5]), find this 
source has a systemic velocity consistent with that of the eclipsing binary. In addition the AO 
observations reveal the 0.7" source is itself a binary. We further discuss the properties of this third 
stellar component in detail in Section [6l 

In order to estimate the brightness of the 3" contaminator and to quantify the contributions 
of additional faint stars 2-6" from KIC 4862625, we use R-band optical imaging from the remote 
SARA Im telescope on Kitt Peak. We obtained a series of 60 60s exposures on 2012 July 18 UT, 
keeping individual exposures short so that all stars of interest remained in the linear count regime. 
The co-added SARA image is presented in Figure [H We used both elliptical model fitting and 
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PSF subtraction (from a similarly bright star without bright neighbors) to remove the wings of 
the bright-star profile. Slight irregularities in the PSF shape made the subtraction of a scaled PSF 
more successful. We estimate errors based on the scatter of photometry in independent fits from 6 
10-exposure averages of the data. The photometric z ero point was establ ished from 20 nearby KIC 



stars brighter than r=17.5, using photometry from iBrown et al.l (|201ll ) .ignoring the small color 
term in transforming from R to Kepler bands. The image quality was 1.5" FWHM, well sampled 
with 0.38" pixels. Measurements used apertures of 4-pixel (1.5"). We estimate an R magnitude of 
18.73 zb 0.05 for the 3" contaminator, providing 1% of the flux measured in the Kepler aperture. We 
also identify an additional faint source located northeast of KIC 4862625. This star is 20.65 it 0.09 
R magnitude, much fainter than the 3" contaminator and contributes well less than a percent flux 
dilution in the Kepler light curve. Any additional sources present within the Kepler photometric 
aperture are fainter than 21st magnitude and will have negligible impact on the Kepler light curve 
for KIC 4862625. 

The dominant sources of additional flux to KIC 4862625 's light curve are the 0.7" and 3" 
contaminators. Combining the contributions from both stars, we assume 12-13% flux contamination 
to the Kepler photometric aperture. Accounting for the dilution of the transit, the flux drop caused 
by the planet now corrected to be 0.1%. We account for this additional light when modeling the 
light curve and assessing PHl's properties (discussed in Section [7]). In Section [9^31 we revise this 
estimate assuming the stellar properties obtained from the photometric-dynamical model and stellar 
evolution models. 



5. Radial Velocity Observations 



Radial veloc ities were obtained with the High Resolution Echelle Spectrometer (HIRES) 
(jVogt et al.l 1 19941 ) on the Keck I telescope on Mauna Kea, Hawaii from 2012 April 1- August 
4. We note that the time baseline and precision of these radial velocity measurements are insuf- 
ficient to constrain the mass of PHI and only characterize the properties of the close binary (Aa 
and Ab). Thes e observations were taken with the stan dard setup used by the California Planet 
Survey (CPS) dHoward et al.1 boiol : Ichubak et all I2OI2I ) with a typical resolution of R ~55,000 
covering 3642 to 7990 A. All Keck observations were taken with the C2 decker, except for the first 
observation, which was taken with the B5 decker. The ob servations were reduced using a standard 



2010 



Chubak et al.ll2012l ;i 



reduction and extraction pipeline developed for the CPS (IHoward et al, 
and then corrected for cosmic rays. 

The radial velocities were measured using the spectral line broadening function (BF) technique 
(jRucinskilll992l ). This method lends itself particularly well to double lined spectroscopic binaries 
where the radial velocity difference between the stars is on the order of the resolution of the spectral 
observation or for high rotation stars where there is significant broadening of the spectral lines. 
To compute the broadening functions, a spectrum of a slowly rotating single star is required as a 
reference template. A single template spectra was used in the BF analysis for all the radial velocity 
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observations of KIC 4862625. The HIRES spectrum of HD 169830, an F7V star, taken on 2012 
April 1 UT was used for this purpose in the BF analysis. The radial velocity of the template star 
was assumed to be -17.4 km s~^. The analysis returns a BF defined over a velocity range for each 
set of spectral lines present in the observation. For each component, the respective broadening 
function is fit to a aanalytic broadening kernel and the peak value is taken as the measured radial 
velocity. A barycentric correction is applied to each of the obtained velocities and the contribution 
of the template radial velocity is removed to obtain the true Doppler shifts for each component 
identified in the BF analysis. 

HIRES spectra are split on three CCD chips: blue (3642.9-4794.5 A), middle (4976.6-6421.4 A), 
and red (6543.2-7989.8 A). The BF analysis is applied separately to the spectra obtained on each 
chip. The majority of HIRES observations were taken with the iodine cell in, which introduces the 
iodine spectra superimposed on the stellar spectra on the middle chip. We find the HIRES spectra 
are double-lined with two BFs: one with a narrow width and a second much broader function. Due 
to the presence of the second velocity component, we chose to restrict ourselves to the blue chip 
to measure radial velocities because the iodine lines made it more difficult to fit the BFs on the 
middle chip, and the BFs on the red chip were weaker than those obtained on the blue chip. A 
representative sample of the BFs obtained from the blue chip is plotted in Figure [71 

Table[3]lists the measured radial velocities and uncertainties for both BFs, and FigureEh plots 
the radial velocities obtained from both the wide and narrow BFs. A two component model was 
used for all observations except for observations 9 and 12, which due to a nearly full moon, required 
an additional third component in order to account for the contribution from sky line emission. The 
peak velocity of the wide component varies as a function of the binary orbital period and is the 
contribution from KIC 4862625. The narrow component is due to another source and is unchanging 
with a nearly constant radial velocity over the 5 month time baseline of the HIRES observations. 
The C2 decker is 14.0" long and 0.861" wide. The B5 decker is the same width but only 3.5" long. 
For all observations the HIRES slit was aligned East- West on the sky. Located at a position angle 
of 123 degrees East of North, the 0.7" contaminator contributes light into the HIRES slit at the 
~10% level, and thus is the additional narrow component identified in the BF analysis. 



The radial velocities were fit with the Eclipsing Light Curve (ELC) code (jOrosz &: Hauschildt 



20001 ). Although the spectra are double-lined the close binary (Aa and Ab) is a single- lined spec- 
troscopic binary (SBl). Star Ab is too faint to be detected in the radial velocity observations, 
and therefore we do not directly obtain the independent masses of Aa and Ab from the HIRES 
observations alone. The ELC modeling ignores the gravitational effects of the planet and assumes 
a two-body system. The ELC modeling cannot constrain the mass of the binary. Instead with the 
radial velocity observations, the ELC fits for the binary's (Aa and Ab) eccentricity (e), inclina- 
tion (i), orbital period (P), time of conjunction (Tconj)) argument of periastron {w), radial velocity 
semi-amplitude (Ka)- The Kepler light curve was used in the fit to force a consistent binary orbital 
phase and eccentricity. The best-fit parameters are reported in Table [H The resulting radial veloc- 
ity model and radial velocity residuals are plotted in Figure [8^. The radial velocity measurement 
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enor for KIC 4862625 was scaled by a factor of 11.2 to produce a reduced of 1 for the model 
fit. We note the 0.7" contaminator has a mean radial velocity (19.12 it 0.49 km s^^) very close to 
that of the systemic velocity obtained for the eclipsing binary, 17.82 it 0.03 km s^^. 

The ELC fit for KIC 4862625 is poor with 1-2 km residuals for the best-fit solution. With 
the CPS observing setup, HIRES is aligned such that that ThAr lines across the echelle format 
fall on the same CCD pixels, to sub-pixel precision such that the obtained HIRES spectra have 
the same wavelength solution, accurate to about a pixel (1.4 km s~^). The HIRES wavelength 
solution zeropoint and dispersion varies within a 1 pixel between observations within a single night 
and on separate nights, which is on the order of 1 km s~^. We do not correct for this effect in 
the BF analysis and this is thus the largest source of uncertainty in our velocity measurements. 
Subtracting the mean radial velocity of the 0.7" contaminator from its measured values, we find 
that the residuals are directly correlated with the ELC fit residuals for KIC 4862625 (see Figure 
[9]). If we assume the 0.7" contaminator has a constant radial velocity, we can use it to correct for 
the instrumental systematics in each observation's wavelength solution. We apply the narrow BF 
velocity residual as a linear offset to correct each of the radial velocities measured for KIC 4862625. 

With the offset-corrected radial velocity observations, the ELC fit is significantly improved 
with typical residuals of 0.5 km s~^ or smaller. The radial velocity measurement error for KIC 
4862625 was scaled by only a factor of 4.63 to produce a reduced of 1- Table [5] lists the offset 
corrected radial velocities for KIC 4862625, and the model parameters are reported in Table [HThe 
ELC model fit model and residuals for the revised KIC 4862625 radial velocities are plotted in 
Figure [Hb- Overplotted on Figure [8h) are the radial velocities of the narrow BF component, the 0.7" 
contaminator, also corrected by the same instrumental offset applied to the KIC 4862625 values. 

With the correction, the systemic radial velocity of KIC 4862625 is found to be 18.18 it 0.03 
km s'^. The average velocity of the 0.7" contaminator (19.12 it 0.49 km s^^) is within the 2- 
a uncertainty of the binary's systemic radial velocity. Within measurement uncertainty the 0.7" 
contaminator has the same radial velocity as the systemic velocity of KIC 4862625. Even before the 
offset corrections for the radial velocity values, the measured systemic velocity of the binary was 
17.82 ± 0.03 km s"\ within 3-a of the mean radial velocity of the 0.7" contaminator. Giv en that 
veloc ity dispersion of random field stars in the Galaxy is between 20 and 60 km s~^ (e.g. IWielen 



19771 ). the radial velocity of the 0.7" contaminator matches the systemic velocity of KIC 4862625. 



At such large separations (~ 1000 AU), the barycentric motion in such a wide binary would be 
negligible and bound stars would have the same proper motions and systemic radial velocities. This 



criter ion has been used to identify wide binary star systems (e.g. iTokovininI 119971 : iPribulla et al 



20091 ). and we adopt it here. Thus Aa, Ab, and the 0.7" contaminator are common radial velocity 



stars and therefore associated. We note that follow-up observations to measure the common proper 
motion or orbital acceleration of these stars in the future could solidly confirm their association. 
Thus, the third source, the 0.7" contaminator, is bound to the eclipsing binary (Aa and Ab) in a 
hierarchical architecture. We further discuss the implications in Section [6l 
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6. Additional Stellar Companions 

The radial velocity observations strongly indicate that the 0.7" contaminator (located at 
roughly ~1000 AU from KIC 4862625) is likely associated with the eclipsing binary (Aa and Ab) 
in a wide binary. Further inspection of the AO observations reveals that the 0.7" contaminator 
is itself a visible double. The point-spread function (PSF) of the 0.7" contaminator is elongated 
compared to the circular PSF of KIC 4862625 (see Figure [TO]l . The one-dimensional full width-half 
maximum (FWHM) in the vertical direction is almost exactly twice as wide as the primary (KIC 
4862625). This was observed on both nights of AO observations. Thus, the 0.7" contaminator is 
itself a wide visual binary, making the KIC 4862625 system a quadruple star system, and the first 
identified to host a planet. 

The AO observations barely resolve the outer star system, which we will label 'Ba' and 'Bb'. 
Our rough estimate of the AO i^s-band brightness difference of 1.0 it 1.0 mag. Both Ba and Bb 
contribute flux to the Keck spectra with nearly equal brightness. Both stars contribute to the second 
broadening function identified in the Keck spectra and both have the systemic radial velocity of 
Aa and Ab. Thus, Ba and Bb are associated with the close binary (Aa and Ab) and are bound 
gravitationally to the system. We are unable to ascertain an accurate projected separation between 
the two stars, but we estimate that it must be less than 0.04" ( < ~40 AU). In Section [931 using 
stellar evolution models, we estimate the spectral type of Ba and Bb and more accurately estimate 
both their projected separation and their projected distance to KIC 4862625. 

This quadruple stellar system is configured as a close eclipsing binary (Aa and Ab) and a wide 
binary (Ba and Bb) arranged with a projected distance between the two binaries that is larger 
than the physical separation between the ind ividu al stars in the binary components. This is not an 



atypical stellar architecture; iTokovininI (|200ll ) and lCorreia et al.l (|2006l ) identify nearly half of their 
sample of quadruple systems in such a two pair hierarchical structure. The two binary pairs are only 
weakly bound to each other. Interactions with passing stars and molecul ar clouds can dynaraically 



disrupt and ionize wide b inaries with separations on the order of 1000 AU (jHeisler fc TremainelflQSg : 



Jiang Sz Tremaindl2010l ). Kaib et al (2012) find for two star wide binaries, after 10 gigayears yu7o 
of systems with semimajor axes greater than 3000 AU are ionized. Those with smaller semimajor 
axes survive much longer. We expect the results to be similar for two -pair hierarchical systems, 
and thus expect such stellar systems separated by a few thousand AU to survive for 10 gigayears 
or longer. Thus the presence of the additional Ba and Bb does not significantly constrain the age 
of the planetary system. 

With the radial velocity observations we can only ascertain that the visual binary is gravita- 
tionally bound to stars Aa and Ab. With only one epoch of adaptive optics observations we do not 
have a sufficient time baseline, to measure their motion on sky, and accurately determine their orbit 
around Aa and Ab. Future adaptive optics observations may be able to measure a pro per motion 



and constrain the orbits of the visual binary. Orbital integrations with MERCURY (IChambers 



19991 ) of Ba and Bb assuming circular orbits with semimajor axes of 1000 AU around Aa and Ab 
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find that on timescales of gigayears, there is httle effect from the distant binary on the planet's 
orbit and the echpsing binary's (Aa and Ab) orbits. This is consistent with the results of Kaib 
et al (2012) who that only for wide-binary systems with semimajor axes less than 10E3 AU, is 
eccentricity pumping observed in the planet's orbit, thus the effect from Ba and Bb on the planet 
is negligible. 



7. Photometric-Dynamical Model 

To constrain the orbital and physical properties of the stars and planet, we employ a 'photometric- 
dynamic al' model, as has been applied to the previously d[iscover ed Kepler transiting circumbinary 



planets (jPovle et al.l 120111 : 1 Welsh et all 120121 : lOrosz et al.l l2012al lbf) fitting the combined radial ve- 



locity observations and Kepler light curve data. In this model, the two host stars (Aa and Ab) 
and the planet interact according to Newtonian gravity. The effects of the wide binary (Ba and 
Bb) are not observable in the Kepler data, and thus were not included in our dynamical model. 
The initial coordinates of the three bodies (Aa, Ab, and the planet) are specified at a particu- 
lar time (in this case to = 2,454,970 BJD) and their positions at any other time are calculated 
via numerical integration. Positions are calculated at the times of the Kepler observations and 
used as inputs to produce a light curve model. The light curve model presumes spherical planet 
and stars with radial brightness profiles specified according to a quadratic limb darkening law 

(/(^) = /o(l - UAa{l - m) - VAa{l " f^?))- 

In previous applications of the photometric-dynamical paradigm, the eclipses of the stellar 
binary and the transit events due to the planet (s) were modeled self-consistently. The goodness- 
of-fit metric, for the best-fitting model, describing the eclipsing stellar binary data alone (i.e., 
considering just the data at primary and secondary stellar eclipses) was verified to be acceptable 
in isolation from the planetary events in these cases. As we will describe in more detail below, the 
joint photometric-dynamical fit to both stellar and planetary eclipse events (and the subsequent 
parameter posterior calculation via Markov Chain Monte Carlo simulation) indicated a solution 
for the stellar events that was marginally unsatisfactory. In detail, the joint fit yielded parameters 
describing the stellar eclipse data that were 1 — 2cr departures from those same parameters when 
fitting the eclipsing binary data alone. In particular, this bias appears to affect the bulk parameter 
estimates of the stars at the ~2-cr level. This error is barely perceptible (by eye) in the photometric 
residuals and only weakly significant in the resulting parameter posteriors; however, the bias is 
statistically clear. We posit that the failure of the model to adequately describe the data is due 
to either missing mass components (e.g., a fourth body), an overly simplified treatment of the 
light curve noise (e.g., correlated noise may be corrupting information gleaned from the planetary 
transits) and/or deterministic light curve features that could not be clearly modeled in this analysis 
(spot crossing during transits). We argue below, after describing our model and analyses, that the 
latter two options are more likely than invoking a fourth body given the significance of the anomaly. 



The photometric-dynamical model describing all eclipse events has 23 parameters. All pa- 
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rameters have uniform priors (over open intervals unless otherwise specified). Eleven parameters 
describe the Jacobian coordinates of the stellar binary and the binary comprised of the center-of- 
mass of the stellar binary and the planet. These eleven parameters arc osculating Keplerian terms 
(one set for each 'binary') that are slowly varying in time: a period, two linearly-independent pa- 
rameters involving the eccentricities and arguments of periastron, a sky-plane inclination, a time 
of conjunction and longitude of the ascending node. Only one nodal angle need be specified as 
the absolute nodal longitude (relative to North) cannot be determined from the data; we fix the 
longitude of the stellar binary to zero. We encode the eccentricity and argument of periastron of 
the stellar binary as components of a vectorial eccentricity: = [e^ sinw^i, cosw^]. Those for 
the binary including the planet involve a square-root in the amplitudes ([y^sinw;,, y^cos w;,]) to 
simplify the specification of uniform priors in and ojb- In addition to the coordinates, we must 
specify the masses of all three bodies. We do so by specifying the mass of star Aa (times the 
gravitational constant), GMAa, and then ratios of Mau- Q = ^Ab/MAa and qp = M^/MAa - both 
of which have positive priors. We parameterize the density of star Aa - this quantity is directly 
measurable from the light curve when the mass ratio q is constrained by the planetary events. 
Six more parameters describe the eclipse and transit light curve: the radius ratios of stars and 
planet, RAb/RAa-, Rb/^Aa, the contaminating flux in the Kepler bandpass relative to that of star 
Aa, Fx/FAa, the flux of star Ab to Aa, FAb/FAa, and two limb darkening parameters of star Aa, 
UAa and VAa- The limb darkening parameters were restricted to physically plausible values with 
the linear constraints < UAa < Ij UAa + VAa < 1- The limb darkening coefficients of the secondary 
star had a negligible affect of the overall fit and were fixed to UAb = 0.6, VAb = 0.1. The width of 
the assumed normal, uncorrelated noise of the photometric data was also fitted [aLc]- Finally, we 
specify the offset of the radial velocity data and a radial velocity jitter term, added in quadrature 
to the formal radial velocity errors. 

Only data within 1 day of either a stellar eclipse or planetary transit were considered. This 
corresponds to 5302 cadences of data at primary or secondary eclipse and 521 cadences of data 
near planetary transit (some of these cadences appear in both data sets as transits occur in close 
proximity to eclipses, at times). These data were initially corrected for a cubic trend in time 
excluding the eclipse events. After an initial photometric-dynamical fit, this cubic fit was repeated, 
this time including all data after having divided out this best-fitting model for the eclipses. The 
level of the out-of-eclipse flux is determined inline prior to the computation of the likelihood (given 
below). In other words, for a given parameter set p, we divide out the eclipse model then determine 
the mean of the residuals and divide through by this value prior to computing the likelihood. This 
linear parameter is weakly correlated with the remaining parameters and trivially marginalized 
over. 

The likelihood of a given parameter set p was specified as 
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X exp 



(Fy/F^g - 0.12)' 
2(0.02)2 



where AF/"^ is the ith photometric data residual, cjlc is the width parameter describing the 
photometric noise of the long cadence data, AVj is the jth radial velocity residual, aj is the 
uncertainty in the jth. radial velocity measurement and any is the stellar jitter term added in 
quadrature with the aj. The final product in the second line of the above equation enforces a 
Gaussian prior on the contamination of Fx/F^a = 0.12 it 0.02, appropriate given the discussion in 
Section HI 

We first performed a fit to the combined eclipse and planetary transit data (hereafter the 'joint 
solution') and determined credible intervals in the above parameters using a Differential Evolution 
Markov Chain Monte Carlo algorithm (DE-MCMC; ted Braak 2008). In detail, we generated a 
population of 100 Markov chains in an overdispersed region covering the final posterior and evolved 
through 19,000 generations until the chains showed adequate convergence according to the standard 
metrics. These chains were concatenated after removing a portion of 'burn-in.' Results from this 
analysis, in the form of medians and 68% credible intervals, are reported in the first data column 
of Table El 

Investigating the parameters from the joint solution, we see that the contamination term is 2a 
lower than the expectation (even with a Gaussian prior enforced). It is not expected that the stellar 
binary eclipse data alone would be sufficient to measure or constrain the level of contamination - 
this drew further investigation. We isolated the data involving the stellar eclipses alone (hereafter 
'EB-only') and repeated our analysis (with the same model). We fixed the planetary mass ratio 
to zero and varied only those parameters relevant to modeling the eclipsing binary (light curve 
and RV data). In particular, we fixed all osculating Keplerian elements of the planetary binary 
and the mass ratio of the stellar binary to those found from the joint solution. Fifty chains were 
evolved through 25,000 generations in the DE-MCMC. The results of this EB-only analysis are 
provided in the second data column of Table El Without any information on the stellar mass ratio 
from the timings of the planetary events, we were unable to provide absolute dimensions to the 
stars or planet; instead, we report standard measured quantities from the analysis of an eclipsing, 
single-lined spectroscopic binary (SBl). This analysis shows, as was expected, that the EB-only 
analysis has no preference in contamination - its marginalized posterior reflects the prior. The 
remaining EB-only parameters, which are correlated with this contamination level, are also 1 — 2o" 
biased relative to the joint solution equivalent parameters as a consequence. The EB-only solution 
was only weakly preferred over the joint solution with A^' ~ 30 when considering just the data 
involving the stellar eclipses, the radial velocity data and the contamination prior (see below for 
further details). 

We suspected that our inability to accurately model the planetary events was driving the 
bias in contamination. To investigate this, we completed a third analysis by fixing the best-fitting 
solution for the eclipsing stellar SBl (from the EB-only analysis above) and fitting for the planetary 
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binary orbit and the stellar mass ratio considering only data involving the planetary events. In 
practice, this involved fixing the 11 osculating elements of the stellar binary, the normalized lengths 
RAa/o.-, RAb/RAai Aux ratios, limb darkening coefficients and the mass function as determined from 
the best-fitting EB-only model. Forty chains were evolved through 64,000 generations. The results 
of this analysis are reported in the final data column of Table [6l In addition to those parameters 
directly measured from this analysis (the osculating elements of the planetary binary, the stellar 
mass ratio and the planetary radius ratio) we report parameter posteriors (daggered values in 
Table ED assuming the independence of EB-only and planetary, EB-fixed solutions. For example, 
the mass of star Aa is determined assuming the mass function determined in the EB-only analysis 
and the mass ratio in the planetary, EB-fixed analysis are independent (simple propagation of errors 
renders our resulting hybrid posterior). 

We compared the best-fitting solutions from the three analyses to help elucidate the nature 
of the anomaly. First, we adopted the value of (Tlc = 203.157 ppm as derived from the EB-only 
analysis. We assert that the model describing the stellar eclipses is less sensitive to potential 
correlated noise sources (stochastic sources or from spot-crossings, for example), as the eclipse 
events are approximately identical at each epoch (outside of insignificant timing anomalies due to 
the planet), and provides a more robust measure of the uncorrelated noise component. We note that 
this value of glc is within 12% of the root-mean-square of the best-fit EB-only residuals. Having 
adopted this width, we compute the traditional statistic for the various contributions to the 
likelihood. The difference in between the joint solution and the EB-only solution considering 
only the radial velocity is negligible. For the 5302 cadences involving only the stellar eclipses 
we find Xe_b ~ 5314.36 and X^^B ~ 5290.02 for the joint and EB-only (or EB-fixed, planetary) 
solutions, respectively. The eclipsing binary solution has 5290 degrees of freedom (after having 
excluded parameters relating to the RV data, the planet, the planet orbit and stellar mass ratio) 
which suggests the EB-only solution is marginally preferred over that from the joint solution. When 
considering the 521 cadences associated with planetary transits alone, we find x^pi = 478.367 and 
XpL ~ 520.572, for the joint and EB-fixed, planetary solutions, respectively. The model for the 
planetary events alone has 512 degrees of freedom. Despite the preference for the lower 'X^p^ from the 
joint solution (which accounts for the overall lower of the joint solution despite the contamination 
prior and the worse-fitting EB solution), the of the EB-fixed solution is acceptable (with a 38% 
probability of having a Xes larger than 520.572). Given this and our confidence in our estimation 
of the contamination prior, we conclude that the preferred EB-only solutions, tabulated in the final 
two columns of Table [6] provide the most accurate assessment of the system parameters. We adopt, 
however, the limit on the mass of the planet from the joint analysis. 



8. ELC Model Fit 



The ELC modeling can be used as an independent check on the binary properties derived 
from the photometric-dynamical modeling effort. The ELC ignores the gravitational effects of the 
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planet, which in this case are neghgible over the duration of the Kepler hght curves. Because Ab 
is too faint to be observed in the HIRES spectrum, ELC cannot obtain an independent mass for 
Aa and Ab or the absolute radii of Aa and Ab. Combining the radial velocity observations and the 
Kepler light curve, the ELC fits for binary's (Aa and Ab) eccentricity (e), inclination (i), orbital 
Period (P), time of conjunction (Tconj); argument of perihelion (w), radial velocity semi-amplitude 
{Ka). the fractional radii of Aa and Ab (R^a/a and R^fe/a), flux contamination, the binary 
temperature ratio (TM/T^Aa), and limb darkening parameters for star Aa (uAa and VAa) assuming 
the same quadratic limb darkening law identical to that used in the photometric-dynamical model 
(as described above in Section [T]). Table [S] lists the physical and orbital parameters for Aa and Ab 
obtained from the ELC fit. The ELC-derived values are consistent with those obtained from the 
photometric-dynamical modeling 'EB only' solution (see Table [6]) within the 1 — a uncertainties. 
Thus we are confident in the physical and orbital parameters derived in the 'fixed EB' solution for 
PHI and its host stars. 



9. Stellar Properties 

In the following section, we further discuss and explore the stellar properties of all four stars 
in the hierarchical system. The resulting best-fitting parameters and their uncertainties for Aa 
and Ab from the photometric-dynamical model are summarized in Table [H Derived physical and 
orbital parameters for Aa and Ab are listed in Table [7] for the fixed EB' solution. The primary star 
(Aa) is an F dwarf ( 1.734 ± 0.044 Rq, 1.528 ± 0.087 Mq) and Ab is M dwarf (0.378 ± 0.023 R©, 
0.408 ± 0.024 M0). Using HIRES spectroscopic observations we further constrain the properties 
of the primary star in the eclipsing binary (Aa). With stellar evolution models, we estimate the 
ages and properties of the planet host stars (Aa and Ab) and the visual binary (Ba and Bb). 



9.1. High Resolution Spectroscopy 



To refine the physical properties of star Aa, the HIRES spectrum from 2012 August 4 UT , 
was fit with Spectroscopy Made Easy (SME) (jValenti Piskunovlll996l : IValenti &: Fisched l2005l ) . 
Modeling the stellar parameters with specific gravity remaining fixed at the value obtained from 
the photometric-dynamical model, log(g)=4.14, the SME analysis finds: Te//=6407 it 150 K, 
T;sim=32.6 ± 2.0 km s-^ and [Fe/H]= +0.21 ± 0.08 for Aa. Although Ba and Bb wih have a 
~10% fiux contribution to the spectrum, the high rotational velocity, which produces a significant 
broadening of the spectral lines, will be the largest source of uncertainty and constraint on the 
precision of the best-fit SME stellar parameters. The temperature estimate is consistent with that 
of an F star. Combining the SME temperature estimate of Aa with the ELC binary temperature 
ratio (Tyift/Tyia), we find Ab has a temperature of 3561 it 150 K. 
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9.2. Stellar evolution modeling 



Combining the results from SME, ELC, and the photometric-dynamical modeling, our de- 
termination of the masses, radii, and temperatures of both stars in the eclipsing binary (Aa and 
Ab), in addition to our knowledge of the metallicity of the system, permit a more stringent com- 
parison with models of stellar evolu t ion with fewer fre e parameters than if the composition were 



unknown (see, e.g., lAndersenI Il99ll : iTorres et al.l 120101 ). This test is particularly interesting for 



the low-mass secondary. Relatively few systems containing M dwarfs with precise (and accurate) 
determinations are available, and these stars have shown signifi cant disagreements with theo ry in 
the sense that they tend to be larger and cooler than predicted (jRibasI 12003 ; Torres &: Ribas 2002 : 
Lopez-Morales k RibaJboOsI : [Morales et allbood : iTorreslbol 2 ) 



Figure [TT] compares the meas ured properties of the Aa star in the Teg-log g plane with models 



from the (Yonsei-Yale) series ( Yi et al.ll200ll ). An evolutionary track (solid line) is shown for 



the measured mass and metallicity of the star, with the 1-cr uncertainty in the location of the track 
coming from the observational errors being represented by the shaded area (dark gray reflecting 
the mass error only, and light gray including also the uncertainty in the metallicity). Isochrones 
are indicated with dotted lines for ages from 1 to 13 Gyr. There is good simultaneous agreement 
with the mass, radius, and temperature of the Aa within the uncertainties, for an age of roughly 2 
Gyr, according to these models. 



In Figure [12] we present another comparison using models from the Dartmouth series (iDotter et al 



20081 ). in which the physical ingredients (including the equation of state and the non-gray bound- 
ary conditions between the interior and the atmosphere) are more realistic for stars significantly 
below 1 Mq. With the reasonable assumption that the composition is the same for the two objects, 
there is again satisfactory agreement between the models and the observations in the mass-radius 
diagram shown in the top panel, also for Ab. The age we infer based on the primary properties 
is 1-2 Gyr, similar to the previous estimate. The predicted flux ratio between the stars in the 
Kp band is 0.0011, according to these models, which is not far from our measurement based on 
the photometric-dynamical model (Table (H). The mass-temperature diagram in the lower panel of 
Figure [12] indicates that the temperatures of both stars are consistent with predictions, given our 
uncertainties, and again suggest an age near 2 Gyr. 

The results for the Ab go against what has usually been found for other low-mass stars in 
eclipsing binaries. Most previous studies have shown that M dwarfs in such systems are typically 
"inflated" by up to 10% or so, and have their temperatures suppressed by up to about 5%. This 
is believed to be due to magnetic activity, as most of those systems are close binaries in which 
the stars are tidally synchronized and are rotating rapidly. This enhances their activity, which is 
usually manifested by modulations in the light curve produced by spots, X-ray emi s sion, and other 



2001 



Chabrier et al 



phenomena, and leads to the observed effects (e.g., iMullan &: MacDonald 
20071 ). For the PHI system we find, on the other hand, that the size of the Ab is consistent with 
theoretical expectations, and so is its temperature. Unfortunately we have no information on the 
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activity level of the Ab, which is some 700 times fainter than the primary. 

The current precision of the measurements is about 6% in the masses and 3-6% in the radii 
of the stars in PHI, while the temperatures are only good to about 150 K, mostly limited by the 
rotationally broadened spectral lines of the primary. Additional spectroscopic observations and 
further data from Kepler should help to reduce some of these uncertainties, providing for an even 
more stringent test of models. 

9.3. Photometric deconvolution 

An additional estimate of the contamination level in the Kp band is possible using available 
multi-color photometry for the system and the stellar evolution models. The photometric measure- 
ments in the KIC, as well as the near-infrared magnitudes from 2MASS, contain significant light 
from four stars: the two components of the eclipsing binary, and those of the 0.7" contaminator. 
We assume here that the 3" contaminator is angularly separated enough that these magnitudes 
are unaffected. With the additional assumption that the four stars are physically associated (see 
Section [6]) , we carried out a deconvolution of all the magnitude measurements using fluxes from a 
model isochrone, for which we adopted a 2-Gyr Dartmouth model that provides a good fit to the 
properties of the eclipsing pair (Aa and Ab)(Sect. l9.2p . Note that those properties (in particular the 
radii and temperatures of the eclipsing components) depend to some extent on the level of contam- 
ination adopted in the photometric-dynamical modeling of the Kepler light curve, so we proceeded 
by iterations between the photometric deconvolution described below, the photometric-dynamical 
modeling, and the model isochrone fitting. 

To deconvolve the combined light of the four stars (SDSS griz, 2MASS JHKs, and the cus- 
tom D51 magnitude from the KIC) we fixed the masses of components Aa and Ab to their values 
obtained from the photometric-dynamical model, which sets their fiuxes in all passbands from the 
isochrone. We then added the masses of Ba and Bb as free parameters, along with the distance 
modulus and reddening, and performed a grid search on those four variables to find the best match 
between the predicted and measured values of the combined light. Appropriate reddening correc- 
tions were applied to each passband. To constrain the masses of stars Ba and Bb we additionally 
used our measurement of the brightness difference of their combined light relative to Aa and Ab, 
from our AO measurements (AJ = 1.89ib0.04, AKg = 1.67ib0.03). As a further constraint to allow 
a distinction between Ba and Bb we used a rough estimate of the Kg-hand brightness difference 
between those stars from the AO image, of 1.0 it 1.0 mag. 

We obtained a good fit to all photometric measurements to well within twice their uncertainties, 
in most cases, and estimates for the masses of stars Ba and Bb of ~O.99M0 and ~O.51M0, 
corresponding approximately to spectral types G2 and M2. The fitted distance modulus corresponds 
to a distance of about 1500 pc. The combined apparent brightness of Ba and Bb in the Kepler band 
amounts to 10.5 it 2.5% of the light of Aa and Ab. From its measured brightness of R = 18.73 it 0.05 
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the 3" companion contributes an additional 1%. This is consistent with our initial estimates of 
contamination presented in Section d] and assumed in our photometric-dynamical modeling. Our 
distance estimate along with the angular separation of the Ba and Bb binary from the target 
translates to a linear separation of about 1000 AU. Stars Ba and Bb are separated by roughly 60 
AU. 



9.4. Stellar variability 

In addition to the eclipses and transits, the Kepler light curve exhibits low-amplitude variability 
on a timescale of a few days. With KIC 4862625 providing ~90% of the flux in the light curve, 
we attribute this variability to Aa. To further analyze this modulation, we mildly detrended the 
light curve by omitting the eclipses and removing low-frequency power by using a polynomial or a 
wide Savitzy-Golay filter (50 cadences on each side). Outliers greater than 4-a from the mean were 
omitted. A Fourier analysis revealed a strong signal with a period of 2.6397 it 0.0014 days; this 
period is consistent with an autocorrelation analysis. The periodicity is present in all Quarters, 
and was confirmed by phase folding each Quarter. The amplitude is on average 300 it 60 ppm, 
though its amplitude is noticeably diminished in Quarters 10 and 11. We note that no correction 
for contamination has been applied in this amplitude estimate. 

In stellar binaries, spin-orbit alignment is expected to be the first to occur. Next is synchro- 
nization of the rotation of th e two stars w ith the orbital motion. Orbital circularization is predicted 



to be the last to occur (see iMazehl |2008| and references within). If we assume the modulation is 
caused by a long-lived cluster of starspots on the primary star (Aa), the periodicity is then the 
rotation period of the star. Using the radius of the star derived from the photometric-dynamical 
model and the inclination of the binary, a predicted projected rotation velocity of f sini = 33.2 it 
0.9 km s^^ is derived, is consistent with the observed value (fsim=32.6 it 2.0 km s^^) obtained 
from the SME analysis. This agreement lends credence to the hypothesis that starspot modulation 
is the source of the periodicity, that the stellar radius is correct, and that the spin axis of the star 
is mostly aligned with the orbital axis of the binary. It follows that that the stellar spin is very far 
from being tidally locke d into a pse udosynchronous spin state; a period of 15.5 days is expected 



using the formulation of iHuti (|l98ll ). This implies the binary (Aa and Ab) is a relatively young 
system with an age less than approximately 9 Gyrs consistent with the findings from the stellar 
evolution modeling in Section 19. 2[ 

We note that we cannot rule out a 7-Doradus type stellar oscillation as the source of the 
2.64 day modulation. This alternative hypothesis is attractive because of the sinusoidal shape of 
the modulation and long-duration stability of the modulation (in amplitude and phase), but the 
change in amplitude in Quarters 10 and 11 and the excellent agreement with the expected rotation 
period from Vrot sin i make the starspot hypothesis more attractive. The primary eclipse fits do 
not show correlated residuals indicative of a starspot crossing, so the starspot hypothesis cannot 
be unambiguously confirmed as was done in the case of Kepler-47 (Orosz et al. 2012b). 
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10. Planet Properties 

Derived physical and orbital parameters for the transiting planet, PHI, are listed in Table 
[71 The best-fit radius for PHI is 6.18 it 0.17 R®. Without eclipse timing variations on the host 
stars (Aa and Ab), the photometric-dynamical model can only place an upper limit on the mass of 
the planet. At the 99.7% confidence level, planet's mass must be less than 169 (0.531 Jupiter 
masses). With a sub- Jupiter mass and radius, the transiting body is well within the planetary 
regime, making PHI Planet Hunters' first confirmed planet. 

Figure [13] shows schematic diagrams of the stellar and planetary orbits. Assuming Aa is the 
dominant source of insolation, neglecting the small impact of Ab on the surface temperature, we 
estimate the equilibrium temperature for PHI. Assuming a Bond albedo of 0.3 (similar to that 
of Neptune) and emissivity of 0.9, we estimate a minimum and maximum attainable equilibrium 
temperatures for PHI. We find a minimum temperature of 524K and a maximum at 613 K. Thus, 
PHI is too hot to be in the habitable zone. Although PHI is a gas giant planet, even if there is 
a possibility of rocky moons orbiting the body, their surfaces would be too hot for liquid water to 
exist. 

The individual planet transits and the best-fitting model are plotted in Figure [T31 All seven 
transits of PHI were across star Aa. With a mutual inclination of 2.619 it 0.057°, the projected 
orbits of PHI and Ab could overlap, but the planet would only transit the face of the Ab near 
quadratures, making this a rare event. If there is a chance alignment where PHI does transit 
across Ab, due to the much fainter magnitude of Ab, the total transit signal-to-noise of that event 
would be below the Kepler point-point noise of 203 ppm. Thus transits across Ab would be at too 
a low signal-to-noise to be identified in the light curve. 

The derived upper mass limit of 169 is the mass required to perturb the orbits of the Aa 
and Ab such that it would induce eclipse timing variations, and therefore likely to be much higher 
than t he tr ue mass of the planet . Using the radius-mass relationship derived from iLissauer et al.l 



(|201ll ) and (iKane fc: Gelindl2012l . private communication), we find a more reasonable mass range 



for PH I lies is somewhere around 20-50 or 0.08 - 0.14 Jupiter masses. As noted by lOrosz et al 



(j2012bl ). all tran siting circumbinary plan ets, including PHI, are Saturn-sized or smaller. This trend 



may support the lPierens &: NelsonI (j2008l ) prediction that Jupiter-mass circumbinary planets should 
be rare due to being ejected or scattered out to the edge of circumstellar disk, but this result is 
far from confirmed. The detection biases have yet to be characterized for any of the circumbinary 
searches, as is the case for sPlanet Hunters as well. PHI was found serendipitously, and a systematic 
search of the Kepler eclipsing binaries has yet to be performed. 
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11. Orbital Stability and Formation 



Following the analytic approximation given bv lHolman Sz WiegertI (|l999l ). the critical orbital 
period and radius in this binary below which the planet's orbit could experience an instability is 
respectively 107 days and 0.55 AU. Like the previously discovered Kepler circumbinary planets, 
with a best-fit orbital period of 138.317 days and semimajor axis of 0.652 AU, the orbit of PHI is 
very close to the edge of stability, being only 29% above the critical perio d and only 18.6% larger 



than the critical orbital radius. Orbital integrations using MERCURY ([Chambersl Il999l ) of the 
system using best-fit physical and orbital parameters for Aa, Ab, and PHI (neglecting the effects 
of Ba and Bb), confirm that the system is indeed stable on gigayear time scales. 

With the long-term stability of the PHI system, the existence of the planet has implications 
for formation scenarios of the stellar system and other multistellar systems. Fragmentation and 
n-body dynamical interactions have been proposed as mechanisms for the formation of stellar 
systems composed of 3 or more stars. The existence of PHI itself suggests that there has been little 
interaction of the outer binary (Ba and Bb) with the inner eclipsing binary (Aa and Ab). Any 
dynamical close encounters between (Aa and Ab) and an additional star early on in its formation 
would have perturbed the surrounding circumbinary disk, likely obliterating the planetesimal disk 
and preventing accreting planetary cores from forming. Thus Ba and Bb have likely not had a 
close encounter with Aa and Ab. The existence of PHI suggests rather than N-body interactions, 
fragmentation in the molecular cloud is the likely scenario for the formation of the four-star system, 
preserving the planetesimal disk around Aa and Ab to grow PHI. 



12. Conclusions 

In this paper we have characterized the properties of the first confirmed planet in the Kepler 
field found by Planet Hunters citizen science project. PHI is a new circumbinary planet orbiting an 
eclipsing binary (KIC 4862625) in a hierarchical quadruple star system. From follow-up observations 
and analysis of the planetary system we conclude the following: 

1. PHI is a 6.18 ± 0.17 R© planet on a 138.317lg;^^7 day orbit. At the 99.7% confidence level, 
we place an upper limit mass of 169 R© (0.531 Jupiter masses). 

2. PHl's host stars are an eclipsing binary (Aa and Ab) with a 20.0002468 ± 0.0000044day 
period. Aa is a 1.734± 0.044 R© and 1.528 ± 0.087 M© F dwarf. Ab is a 0.378 ± 0.023 R© 
and 0.408 ± 0.024 M© M dwarf. 

3. Adaptive optic observations reveal an additional pair of stars in a visual binary (Ba and Bb) 
contaminating the Kepler photometric aperture, likely composed of a G2 and a M2 stars, at a 
projected distance of simlOOO AU. Radial velocity observations confirm the pair is associated 
and bound to the planetary system. 
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4. We estimate the age of the planetary system based on stellar evolution models of Aa and Ab 
at ~2 Gyrs. 



The frequency and range of orbital configurations for circumbinary and other multistar plane- 
tary systems has yet to be fully explored, but with this seventh confirmed transiting circumbinary 
planets, we are moving closer to probing this parameter space. These planetary systems provide 
new boundary conditions for planet formation This ever increasing sample of dynamically extreme 
planetary systems will serve a s unique and vital tests of proposed planetesimal formation mod- 



els and core accretion theories ( Bate et al. 



200C 



Pierens &: Nelson 



2008 



Paardekooper et al. 



2012 



Meschiaril 12012 : Ulexandeil boid ). PHI was found serendipitously by Planet Hunters volunteers 
examining the light curves of the known Kepler eclipsing binaries. The discovery of PHI highlights 
the potential of visual inspection through crowd sourcing to identify planet transits in eclipsing 
multistar systems, where the primary and secondary eclipses of host stars make detection challeng- 
ing compared to planets with single host stars. In particular, Planet Hunters may exploit a niche 
identifying circumbinary planetary systems where the planet is not sufficiently massive to produce 
measurable eclipse timing variations that would be identifiable via automatic detection algorithms. 
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Table 2:: Measured primary and secondary eclipse times for 
KIC 4862625 



Orbital 


Primary Eclipse 


Secondary Eclipse 


Cycle 


Midpoint 


Uncertainty 


Midpoint 


Uncertainty 




(BJD- 2455000) 


(minutes) 


(BJD- 2455000) 


(minutes) 


1 


-32.18064 


0.889 


-24.32048 


7.736 


2 


-12.18029 


0.889 


-4.31919 


6.602 


3 


7.82051 


0.889 






4 


27.82051 


0.876 


35.67367 


8.872 


5 


47.81987 


0.967 


55.68321 


8.682 


6 


67.82051 


0.903 


75.66013 


7.454 


7 


87.82091 


1.032 


95.68077 


8.493 


8 


107.82066 


0.850 


115.67901 


7.169 


9 


127.82052 


0.863 


135.67771 


5.848 


10 


147.82173 


0.993 







Continued on Next Page. . . 



This preprint was prepared with the AAS IM^^X macros v5.2. 



Table 1: Measured trans it times for PHI 

Quarter Transit Midpoint 1-a Uncertainty 





(BJD- 2455000) 


(Minutes) 


2 


70.80674 


8.94 


4 


207.44916 


9.00 


5 


344.11218 


6.14 


7 


479.98746 


8.89 


8 


613.17869 


9.04 


10 


749.20260 


6.49 


11 


885.91042 


8.94 
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Table 2 - Continued 



Orbital 


Primary Eclipse 


Secondary Eclipse 


Cycle 


Midpoint 


Uncertainty 


Midpoint 


Uncertainty 




(BJD- 2455000) 


(minutes) 


(BJD- 2455000) 


(minutes) 


11 


167.82246 


0.850 


175.67458 


6.885 


12 


187.82161 


0.876 


195.67849 


6.791 


13 


207.82279 


1.229 


215.67642 


7.169 


14 


227.82271 


0.889 


235.67813 


6.696 


15 


247.82347 


0.889 


255.67295 


7.358 


16 


267.82365 


0.876 


295.68469 


7.074 


17 


287.82313 


0.889 


315.66953 


7.839 


19 


327.82365 


0.902 


335.67685 


6.885 


20 


347.82402 


0.889 


355.68154 


7.764 


21 


367.82395 


0.902 


375.67816 


7.35804 


22 


387.82389 


0.889 


395.68198 


6.980 


23 


407.82594 


0.902 


415.68345 


6.885 


24 


427.82464 


0.876 


435.68501 


6.980 


25 


447.82520 


0.915 


455.67885 


6.980 


26 


467.82607 


0.889 


475.68108 


7.169 


27 


487.82642 


0.889 


495.67777 


7.641 


28 


507.82580 


0.941 


515.67622 


7.930 


29 


527.82476 


1.006 


535.68698 


7.074 


30 


547.82560 


0.889 


— 


— 


31 


— 


— 


575.68315 


7.454 


32 


587.82759 


0.889 


— 


— 


33 


607.82824 


1.006 


615.67919 


6.696 


34 


627.82612 


0.876 


655.69609 


7.925 


35 


647.82747 


0.876 


675.69147 


7.650 


36 


667.82816 


0.889 


695.68003 


7.263 


37 


687.82937 


0.889 


715.67600 


7.454 


39 


727.82918 


0.993 


735.68862 


7.547 


40 


747.82850 


1.006 


755.68434 


6.980 


41 


767.82945 


0.876 


775.68358 


7.263 


42 


787.83000 


0.967 


795.67857 


8.211 


43 


807.82964 


0.889 


815.68275 


6.980 


44 


827.83067 


0.902 


835.68385 


7.265 


45 


847.82987 


0.902 


855.68746 


6.696 


46 


867.82980 


0.876 


875.68666 


8.114 


47 


887.83118 


0.902 


895.68856 


7.547 
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Table 2 - Continued 



Orbital 


Primary Eclipse 


Secondary Eclipse 


Cycle 


Midpoint Uncertainty 


Midpoint Uncertainty 




(BJD- 2455000) (minutes) 


(BJD- 2455000) (minutes) 


48 


907.83190 0.889 


915.68830 6.696 


49 


927.83150 0.915 
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Table 3: Measured Radial Velocities 









Broad BF 


Narrow BF 




Observation 


UT Date 


BJD-2455000 


Radial 




Radial 


0"RV 


# of 


# 






Velocity 


(km s^^) 


Velocity 


(km s~^) 


BFs 








(km s ^) 




(km s ''") 




Fit 


± 


9019 Anril 1 


-LU -Li7. 1 -LOU / OU 


34.647 


0.111 


18.935 


0.189 


9 


2 


2012 May 26 


1073.9616350 


-0.050 


0.186 


19.299 


0.256 


2 


3 


2012 May 28 


1076.0832280 


19.503 


0.126 


19.963 


0.218 


2 


4 


2012 June 19 


1098.0375080 


32.848 


0.176 


20.736 


0.367 


2 


5 


2012 June 19 


1098.1085180 


32.989 


0.153 


22.249 


0.290 


2 


6 


2012 June 20 


1099.0462140 


31.941 


0.116 


16.303 


0.333 


2 


7 


2012 June 20 


1099.0574410 


31.925 


0.094 


16.722 


0.176 


2 


8 


2012 July 1 


1110.0675670 


2.359 


0.140 


20.088 


0.133 


2 


9 


2012 July 2 


1111.0639150 


-4.437 


0.143 


17.350 


0.262 


3 


10 


2012 July 3 


1112.0969640 


-8.515 


0.153 


16.860 


0.268 


2 


11 


2012 July 4 


1113.0660410 


-5.516 


0.275 


19.300 


0.394 


2 


12 


2012 August 4 


1144.0612020 


30.766 


0.133 


20.612 


0.196 


3 



Table 4: Preliminary ELC Radial Velocity Fit Parameters for KIC 4862625 



Parameter 



Uncorrected Radial Velocities Corrected Radial Velocities 



Ka (km s-i) 
P (days) 

Tconj (BJD- 2455000) 

ecos(u)) 

esm{w) 

e 

w (degrees) 



21.490± 0.480 
20.0002514 ± 0.0000056 
-32.18068 ± 0.00010 
-0.16739 ± 0.00006 
-0.14622± 0.00093 
0.2223 ± 0.0010 
221. 14± 0.19 



21.185 ± 0.205 
20.0002489 ± 0.0000051 
-32.18068 ± 0.00099 
-0.16763 ± 0.00010 
-0.13363± 0.00482 
0.2144 ± 0.0029 
218.56 ± 1.04 
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Table 5: Corrected Radial Velocities for KIC 4862625 



Observation 


UT Date 


BJD-2455000 


Radial Velocity 


0"RV 


# 






(km s^^) 


(km s~^) 


1 


2012 April 1 


1019.1136730 


34.746 


0.111 


2 


2012 May 26 


1073.9616350 


-0.314 


0.186 


3 


2012 May 28 


1076.0832280 


18.575 


0.126 


4 


2012 June 19 


1098.0375080 


31.146 


0.176 


5 


2012 June 19 


1098.1085180 


29.775 


0.153 


6 


2012 June 20 


1099.0462140 


34.672 


0.116 


7 


2012 June 20 


1099.0574410 


34.237 


0.094 


8 


2012 July 1 


1110.0675670 


1.306 


0.140 


9 


2012 July 2 


1111.0639150 


-2.752 


0.143 


10 


2012 July 3 


1112.0969640 


-6.340 


0.153 


11 


2012 July 4 


1113.0660410 


-5.782 


0.275 


12 


2012 August 4 


1144.0612020 


29.189 


0.132 



Table 6:: Model parameters for PHI. The reference epoch 
is to =2,454,970 (BJD). f - parameter median and credible 
intervals determined from a combination of EB-only and EB- 
fixed, planetary solutions. 



PariUU(;Ler Nmuk; 



Joint S(jluliou 



EB oulv 



Plaucn. only (EB Fixed) 



Mass parameters 

Mass of Star Aa, M^a (Mq) 

Mass Ratio, MAb/MAa 

Mass of Star Ab, M^b (Mq) 

Planet Mass Ratio, M^/MAa (99.7% conf.) 

Planet Mass, Mb (M®, 99.7% conf.) 
Radius Parameters 

Normalized Radius, RAa/o,A 

Radius of Star Aa, Rau (Rq) 

Radius Ratio, RAb/RAa 

Radius of Star Ab, RAb {Rq) 

Planet Radius Ratio, Rb/RAa 

Planet Radius, Ri, (i?©) 
Planetary Orbit 

Period, Pb (days) 

Eccentricity Parameter, ^yebCos{u!b) 
Eccentricity Parameter, y^sin(a;b) 
Bary. Transit Time, tb (BJD-2,454,900) 
Sky-plane Inclination, (deg) 
Relative Nodal Angle, (deg) 
Stellar Orbit 
Period, Pa (days) 
Eccentricity Parameter, 6^005(0;^) 
Eccentricity Parameter, ca sm{ijOA) 



1.384 ± 0.079 
0.2794 ±0.0051 
0.386 ± 0.018 

< 0.000357 

< 169. 

0.04690 ± 0.00041 
1.759 ±0.030 



0.240 



h0.020 
-0.015 

422+°-°3^ 

U.^ZZ_Q Q2g 

0.03222 ± 0.00045 
6.18 ±0.14 



138.506l°:J°^ 
0.228lS;°iI 

-0.044lHi 
174.71 ±0.11 
yu.uzz_Q 
0.89 ±0.16 



on non2i4+°-°°°°^-^ 

ZU.UUUZi^_Q QQQQ43 

n 1 «77q+0.00021 
U.ID/ 'O_o.00018 

-0.1291 ± 0.0084 



(fixed) 

0.04490 ± 0.00064 
0.219 ±0.014 



20.0002468 ± 0.0000044 
-0.16747 ± 0.00037 
-0.141 ±0.015 



1.528 ±0.087t 

U.ZUUO_g QQ23 

0.408 ± 0.024t 
(fixed) 



0.04516 (fixed) 
1.734 ±0.044t 
0.215 (fixed) 
0.378 ± 0.023t 
0.03268 ± 0.00032 
6.18±0.17t 



+0.040 

-0.027 
0.0044 
0.0060 

0.0824^°;°°^! 



138.317 
0.2519+ 



174.490lHi 
90.050 ± 0.053 
0.89 ±0.14 

20.00024678 (fixed) 
-0.16750 (fixed) 
-0.1393 (fixed) 



Continued on Next Page. 



Pariimcler Naiuc; 


Joiul S(jluliou 


EB ouly 


Plaucn. only (EB Fixed) 


Primary Eclipse Time, tA (BJD-2,454,900) 


67.81776 ± 0.00015 


67.81797 ± 0.00016 


67.81798 (fixed) 


Sky-plane inclination, ijy (deg) 


87.360in?2 


87.592 ± 0.053 


87.58 (fixed) 


Light Parameters 








Linear Limb Darkening Parameter, UAa 


0.27 ±0.15 


^•"-'-'--0.21 


0.32 (fixed) 


Quadratic Limb Darkening Parameter, VAa 


0.08 ± 0.18 


n 17+0.21 

-0.28 


0.17 (fixed) 


Stellar Flux Ratio, FAh/FAa xlOO) 


1 t;QQ+0.0142 

U.±000_Q QQgg 


n 1 41 4+0-0040 

U.14H41:_Q QQ3g 


0.1405 (fixed) 


L/ontammation, rx/i^Aa (xiUUj 


7.0 ± i.y 


11.8 ± 1.9 


12. (lixedj 


Radial Velocity (RV) Parameters 








RV Semi-amplitude, KAa (km/s) 


20.69 ±0.31 


20.56 ±0.33 




Mass Function, f{MAa,MAb) {Mq) 


0.01836 ± 0.00082 


0.01802 ± 0.00089 


0.01807 (fixed) 


RV Jitter (km/s) 




0.8ll°:?3^ 





Table 6:: Model parameters for PHI system. The reference 
epoch is to =2,454,970 (BJD). f - parameter median and 
credible intervals determined from a combination of EB-only 
and EB-fixed, planetary solutions. 



Table 7: Derived planet and binary parameters from photometric-dynamical model fit. Reference 

epoch is to =2454970 (BJD) . 

Parameter Name Joint Solution EB only Planet only (EE Fixed) 

Stellar (Aa &: Ab) Bulk parameters 

Density of Star Aa, pAa (g/cm^) 0.2542 ± 0.0076 - 0.293 ± O.OlSf 

Density of Star Ab, pAb (g/cm^) 5.1 ± 1.2 - 

Gravity of Star Aa, log(p^„) (cgs) 4.089 ± 0.014 - 4.144 ± 0.014t 

Gravity of Star Ab, log{gAb) (cgs) 4.772+°;°^^ - 4.892 ± 0.057t 
Stellar (Aa & Ab) Orbital Properties 

Semimajor axis , qa (AU) 0.1744 ± 0.0031 - 0.1797 ± 0.0035t 

Eccentricity, ca 0.2117 ± 0.0051 0.2188 ± 0.0094 0.2179 (fixed) 

Argument of Periastron ua (deg) 217.6 ± 1.9 220.ll:|^ 219.7504 (fixed) 
Planet Orbital Properties 

Semimajor axis Ub (AU) 0.634 ± 0.011 - 0.652 ± 0.012t 

Eccentricity of Planetary Orbit, Cb 0.0539 ± 0.0081 - 0.0702;°;°°^9 

Argument of Periastron, Wfe (deg) 348.0l||^ - 341.88to:?o 

Mutual Inclination", /(deg) 2.814 ± 0.073 - 2.619 ± 0.057 



"Stellar to planetary orbit mutual inclination: The mutual inclination is the angle between the orbital planes of the 
binary and the planet, and is defined as cos / = sin ia sin if, cos Af2 + cos ia cos it. f - parameter median and credible 
intervals determined from a combination of EB-only and EB-fixed, planetary solutions. 



Table 8. 



Parameters from full ELC Model 



Parameter Best fit 

Ka (km s-^) 21.179 ± 0.162 

e 0.2165 ±0.0013 

uj (dcg) 219.30 ± 0.45 

RAa/a 0.04489 ± 0.00007 

RAb/a 0.00938 ± 0.00008 

TeS,Ah/TeSM 0.5557 ± 0.0062 

i (deg) 87.612 ± 0.004 

it^a 0.208 ± 0.020 

w^a 0.293 ± 0.0249 

P (days) 20.0002490 ± 0.0000049 

Tconj (BJD-2455000) -32.18064 ± 0.00203 

Contamination 0.109 ± 0.012 



Note. — Note: Subscript "Aa" denotes the 
primary star, subscript "Ab" the secondary star 
in the eclipsing binary KIC4862625 . 
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Fig. 1. — Planetary transits in addition to binary-star eclipses identified by Planet Hunters volun- 
teers in KIC 4862625 as viewed in the Planet Hunters interface. 
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Fig. 2. — Top: Normalized Kepler light curve for KIC 4862625 for Quarters 1-11 Bottom: Individual 
transits of the planet across star Aa identified in Quarters 1-11. In Quarter 4 (Q4), the planet transit 
occurs just before a primary eclipse. 



- 38 - 




200 400 600 800 

Time (BJD-2455000) 

Fig. 3. — 0-C offsets measured for the primary and secondary (red in the online version) eclipses 
of KIC 4862625. 
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Fig. 4. — 0-C offsets and transit durations for the circumbinary planet, a) 0-C offsets measured 
for tlie planet transits versus time. The l-cr error bars are plotted. With uncertainties on the 
transit times on the order of ~5 minutes, the error bars are smaller than the size of the symbol, 
b) 0-C offsets measured for the planet transits versus binary orbital phase. The l-cr error bars 
are plotted. With uncertainties on the transit times on the order of ~5 minutes, the error bars 
are smaller than the size of the symbol, c) Transit duration as a function of binary orbital phase. 
Transits near secondary eclipse (p = 0-39 are wider than those close to primary eclipse (j) = 0.0 
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Fig. 5.— NGS AO Ks image of KIC 4862625 and the 0.7" contaminator located 0.702" away from 
the primary Kepler target. 
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Fig. 6. — SARA 0.9m stacked image of KIC 4862625 and surrounding stars. The inner circle (light 
blue in the online version) has a radius of 6" and the outer circle (red in the online version) has a 
radius of 10" centered on KIC 4862625. 




Fig. 7. — The Broadening functions and the analytic broadening kernel fits (solid line for the 
sum, dotted lines for each component separately) for a subset of the Keck HIRES radial velocity 
observations taken at different binary orbital phases (i;^). 
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Fig. 8. — a) Top: Uncorrected radial velocities obtained for KIC 4862625 and ELC fit to the light 
curve and radial velocity observations. Bottom: The residuals of the ELC fit. b) Top: Offset- 
corrected radial velocities obtained for KIC 4862625 and improved ELC fit to the light curve and 
radial velocity observations. Bottom: The residuals of the ELC fit. - In both plots: The filled 
dots are the measurements for KIC 4862625. The circles are the radial velocities of thenarrow 
broadening function component, the 0.7" contaminator. Solid line (red in online version) is the 
best-fitting ELC solution. 
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KIC 4862625 residuals (km s"') 

Fig. 9. — Radial velocity residuals for KIC 4862625 versus the radial velocities residuals of the 
narrow broadening function component, the 0.7" contaminator before applying the offset-correction. 
The plotted line (red in the online version) is the 1 to 1 correlation in primary and companion radial 
velocity residuals. 
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Fig. 10. — Contour plot of the 0.7" contaminator (starting at 1049.404 counts) and KIC 4862625 
(starting at 3940.448 counts) from the stacked NGS AO NIRC2 Ks observations. Contour intervals 
are 300 counts in both plots. 
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Fig. 11. — Stellar evolution models from the (Yonsei-Yale) series by lYi et al.l ( 200ll ) compared 
against the measurements for the Aa. The solid line represents an evolutionary track for the 
measured mass (MAa = 1.528 ±0.087 Mq) and metallicity ([Fe/H] = +0.21 ±0.08), and the dotted 
lines are 1-13 Gyr isochrones (left to right) in steps of 1 Gyr. The uncertainty in the location of 
the track that comes from the observational errors is shown by the shaded area. The darker area 
reflects the mass error, and the lighter area includes also the uncertainty in the metallicity. 
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Fig. 12. — Mass, radius, and effective tempera tures for Aa a nd A b compared against 1-13 Gyr 
model isochrones from the Dartmouth series by iDotter et al.l (|2008l ) , for the measured metalhcity 
of the system ([Fe/H] = +0.21 zt 0.08). The solid line corresponds to a 2 Gyr isochrone, for 
reference. (Top) Mass-radius diagram showing good agreement for both stars (Aa and Ab) for an 
age between 1 and 2 Gyr; (Bottom) Mass-temperature diagram again showing that the models 
match our estimated temperatures, within the errors. 
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Fig. 13. — Orbital configuration of the PHI system for the reference epoch, to =2,454,970 
(BJD).The inner close binary (primary star Aa and the secondary star Ab) and the planet (la- 
beled b) are shown. The outer visual binary (Ba and Bb) is not depicted. Top: A scaled, face-on 
view of the orbits. On this scale the stars and the planet are too small to be seen, and their 
positions are represented by the boxes. Bottom: The region between the vertical lines in the top 
diagram is shown on an expanded scale with an orientation corresponding to what would be seen 
from Earth. 
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Fig. 14. — Photometric-dynamical model fit to individual planet transits. In Quarter 4 (Q4) , the 
planet transit occurs just before a primary eclipse. 



